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ABSTRACT 


This thesis explores how an Asteroid Sample Return Mission might make 
use of solar electric propulsion to send a spacecraft on a journey to the asteroid 
1989ML and back. It examines different trajectories that can be used to get an 
asteroid sample return or similar spacecraft to an interplanetary destination and 
back in the most fuel-efficient manner. While current plans call for keeping such a 
spacecraft on the asteroid performing science experiments for approximately 90 
days, it is prudent to inquire how lengthening or shortening this time period may 
affect mission fuel requirements. Using optimal control methods, various mission 
scenarios have been modeled and simulated. The results suggest that the 
amount of time that the spacecraft may spend on the asteroid surface can be 
approximated as a linear function of the available fuel mass. Furthermore, It can 
be shown that as maximum available thrust is decreased, the radial component 
of the optimal thrust vector becomes more pronounced. 
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I. INTRODUCTION 


A. THE HISTORICAL PERSPECTIVE 

It is frequently asserted that the history of modern space propulsion began 
with the invention of gunpowder in China over a thousand years ago. It is around 
this time that the first solid fuel chemical rockets appeared. Even with the later 
development of liquid fueled rockets and space launch vehicles in the Twentieth 
Century, many of the fundamental principles that govern the construction and 
operation of modern rocket propulsion systems have changed little since the time 
when ancient Chinese rocket engineers like the legendary Wan Hoo are said to 
have attempted the first space flight using a rocket propelled wicker chair. Those 
present for the occasion, and prudent enough to have stood well clear of the 
would be astronaut, were among the first to fully understand and appreciate the 
tremendous unwieldy power and volatility of the rocket, to say nothing of the 
difficulties inherent in applying any form of control to the process of thrusting 
[Ref. 14]. A similar understanding of rockets combined with healthy respect for 
their limitations has persisted down the centuries and continues to inform every 
aspect of space propulsion engineering to this day. 

In spite of the risks, rockets continue to be the preferred method of 
primary propulsion. Understandable, since a considerable part of any modern 
day space mission involves simply getting into orbit. Rocket propulsion works in 
the Earth’s atmosphere more or less as well as it does in space. Moreover, it can 
be tested and developed here on Earth without the use of expensive vacuum 
chambers. Given the fledgling state of current space fairing capability, chemical 
rocket propulsion makes perfect sense. 

In the near future however, space mission planning is increasingly likely to 
begin not at launch, but with a highly maneuverable spacecraft already placed in 
orbit. Mission planning, to the extent that it occurs prior to launch, is likely to be 
open ended with specific trajectories left undefined until the need arises. In many 
instances, the inherent risks and inefficiencies of the chemical rocket will make it 
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an unattractive alternative for primary propulsion. Space propulsion systems of 
the future will combine unprecedented efficiency with a degree of reliability, 
flexibility, and control not seen in modern chemical rockets. Designed to function 
only in the vacuum of space, their impact on exploration is likely to be as 
profound as that of the innovations in sailing technology that preceded the 
Sixteenth Century voyages of discovery. 

B. THE EMERGENCE OF ION PROPULSION 

In October of 1998, a spacecraft called Deep Space One was launched 
from Cape Canaveral. Its primary mission was to test a number of high-risk 
technologies. Foremost among these was the Boeing NSTAR 30 cm gridded ion 
thruster. Though ion thrusters had been around for decades. Deep Space One 
was the first spacecraft to successfully employ an ion thruster for primary 
propulsion on an extended deep space mission [Ref. 12]. 



Figure 1. Deep Space One [From: Ref. 1] 


Over the next three years, NSTAR would perform flawlessly over 16,000 
hours of continuous operation [Ref. Ij. Since then, the European Space Agency 
has also successfully employed an ion thruster to propel Europe’s first lunar 
mission, known as SMART 1, into orbit around the moon. Flaving proven itself to 

2 




be a safe, effective, and reliable engine for primary propulsion, the ion thruster is 
presently being integrated into many spacecraft now under development. 

Ion thrusters work by using electricity to ionize an inert gas such as xenon 
and accelerate it out of the spacecraft at extremely high velocity. While actual 
thrust is imperceptibly low when compared to chemical rockets, ion thrusters are 
capable of applying this thrust continuously for months at a time. Fuel mass 
efficiency, as measured in terms of specific impulse, can be more than ten times 
higher than that of chemical rockets. 

While a spacecraft thus propelled can frequently take much longer to 
reach a particular destination than a conventional chemical rocket, it need only 
carry a small fraction of the chemical thrusters’ required fuel mass, thus leaving 
more room for additional payload. Moreover, an ion thruster employs no moving 
parts and uses a chemically inert propellant. Energy is not stored in the form of 
volatile chemical explosives but is extracted from the environment through the 
spacecraft’s solar panels. This approach gives the ion thruster a tremendous 
safety advantage over chemical thrusters. 

C. A PROBLEM WELL SUITED TO OPTIMAL CONTROL 

In addition to these advantages, ion thrusters are much easier to control 
than chemical thrusters. On a fundamental level a chemical thruster propels a 
spacecraft by applying a series of prolonged explosions. The effect on spacecraft 
trajectory cannot always be precisely predicted and this uncertainty must be 
accounted for in the planning of spacecraft maneuvers. Uneven burning of solid 
propellants frequently necessitates spin stabilization of the entire spacecraft 
while a burn is underway. Throttling or shutting down a solid booster before all 
the fuel is exhausted is impracticable. Liquid chemical thrusters can be throttled, 
shutdown, and restarted. But the repeated cycling of valves over time can 
contribute to increasingly uncertain performance in the long term. 

Ion thrusters, on the other hand, can be throttled up and down to the 
nearest thousandth of a Newton merely by changing the applied electrical power. 
While thruster grid erosion can result in slight performance degradation over 
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time, the absence of moving parts enhances reliability and ensures predictable 
performance over the life of the spacecraft. Lengthy burn times equate to long 
reaction and recovery times that facilitate timely intervention if thrusting is not 
having the desired effect on spacecraft trajectory. Moreover, thruster pointing is 
easily accomplished with conventional momentum wheels or gimbaled mountings 
even while thrusting is at maximum power. 

Such enhanced control and flexibility facilitates a wide variety of useful 
tasks including precise maneuvering for rendezvous and surveillance of asteroids 
or other spacecraft. It makes docking for servicing, refueling, and repair much 
safer and easier to control and automate. It also makes possible the kind of 
precise formation flying necessary for sparse aperture remote sensing using 
multiple spacecraft. 

Enhanced controllability naturally translates into greater opportunity to 
exercise optimal control. Spacecraft trajectories can be designed to minimize 
transit time, minimize fuel consumption, maximize precision, or optimize some 
other desired mission parameter. Trajectories can even be designed to optimize 
a precisely defined combination of parameters in a precisely defined way - a 
problem that naturally lends itself to real time computerized automation. 

D. THE FUTURE OF SPACE PROPULSION 

Not surprisingly, government organizations and private companies from 
around the world are even now developing and marketing ion propulsion 
spacecraft to fill a wide variety of roles and missions. In the UK, the recently 
privatized consortium of Defense Ministry laboratories known as Qinetiq PLC is 
marketing a line of highly maneuverable microsatellites to provide space 
situational awareness to its customers including the ability to rendezvous with 
and photograph other spacecraft in geosynchronous orbit. The spacecraft weighs 
less than two hundred kilograms and uses a single ion thruster for primary 
propulsion [Ref. 2]. 

Also in Europe, Orbital Recovery Corporation is marketing the 
ConeXpress Orbital Life Extension Vehicle (CX OLEV) as an ion propulsion 
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“space tug” to facilitate the refueling and servicing of communications satellites 
that have expended their onboard propellant supplies but are otherwise quite 
functional, thereby giving a new lease on life to a very expensive capital 
investment. The company foresees a large market among the global 
telecommunications industry. In a stroke of innovative engineering, the 
spacecraft structure actually doubles as the Ariane V payload adapter; thus, 
preserving valuable payload volume for paying customers [Ref. 3]. 



Figure 2. ConeXpress Orbital Life Extension Vehicle [From: Ref. 3] 



Figure 3. Jupiter Icy Moons Orbiter [From: Ref. 4] 
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Large arrays of ion thrusters are envisioned as the primary propulsion 
system for NASA’s Jupiter Icy Moons Mission. Power is to be provided by an 
onboard nuclear reactor under development as part of NASA’s Project 
Prometheus [Ref. 4]. 

In the near term, NASA is proposing a number of smaller deep space 
missions to be propelled by ion thrusters. An asteroid sample return mission, the 
subject of this thesis, is one of these. But there should be no doubt that the utility 
and applicability of the ideas discussed here in the context of the asteroid sample 
return mission extend well beyond any single spacecraft or any single type of 
mission. Whether for military, commercial, or purely scientific purposes, optimal 
control of low thrust trajectories is likely to become an increasingly important 
aspect of space mission planning in the years ahead. 
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II. BACKGROUND 


A. THE ASTEROID SAMPLE RETURN MISSION 

The Near Earth Asteroid Sample Return Mission is intended to provide 
planetary geologists an opportunity to examine asteroid mineral samples up 
close for the first time. In February 2001, the NEAR spacecraft successfully 
landed on the asteroid Eros and provided detailed data on chemical and mineral 
composition as well as numerous close up photographs of the asteroid’s surface 
features. While this kind of information can be of immense value in answering 
fundamental questions about the nature of the solar system, actual samples of 
asteroid surface material could yield still greater discoveries by allowing 
scientists to test the material back on Earth [Ref. 5]. 


I 5 mdu Fixed 
High Gam Antenna 



Solar Array G 





n’t/ Batteries 


Solar Array Release 
Mechanism 



Figure 4. Reference Sample Return Spacecraft [From Ref. 5] 
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B. TARGET AND TRAJECTORY 

Asteroid 1989ML is the reference mission destination for this thesis. But 
among the hundreds of possibilities, there are many attractive candidates. So 
much so that a multiple asteroid rendezvous has been considered. Deciding on 
which bodies to visit and in what sequence is itself a difficult optimization problem 
[Ref. 17]. 

Figure 5 below shows two circles approximating the orbit of the Earth 
(blue) and Mars (red) based on their semi-major axes. 1989ML is roughly one 
kilometer across and lies in a slightly eccentric orbit (green) between the two. 
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C. THE ASTEROID SAMPLE RETURN SPACECRAFT 

The requirement to actually fly out to the asteroid and return poses a 
number of technical challenges beyond those encountered in a one way journey 
such as that undertaken by NEAR in 2001. Spacecraft mass becomes an even 
more critical issue since the spacecraft must now carry sufficient fuel for a round 
trip. The mission can be expected to last at least twice as long, thus mandating 
high component reliability and reduced risk, qualities not frequently associated 
with pressurized tankage, valves, and piping systems. The target asteroid, 
1989ML, being less than one kilometer in diameter, exerts very little gravitational 
pull. Precise maneuvering in such a weak gravity field necessitates a very fine 
level of thrust control. The ion thruster favorably addresses all of these issues 
and is thus particularly well suited to fulfill the primary propulsion requirements 
for this mission. The reference mission uses three NSTAR ion thrusters; two for 
propulsion and one spare. Each thruster will be similar to the one flown for the 
first time on Deep Space One. 
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III. FORMULATING THE ASTEROID SAMPLE RETURN 

PROBLEM 


A. MODELING THE DYNAMICS 

1. Two Dimensional Model 

Since the orbital inclination of 1989ML is only 4.38 degrees, the sample 
return mission can be modeled two dimensionally using the small angle 
approximation. This greatly simplifies the calculations while still yielding 
reasonably accurate results. [Ref. 7] 

2. Coordinate System, State Vector, and Control Vector 

The two dimensional model naturally lends itself to a heliocentric polar 
coordinate system in which the spacecraft state vector can be described in terms 
of radial position r, angular displacement with respect to Aries 6, radial velocity 
u , transverse velocity v, and spacecraft mass m : 


r 


e 


X - 


u 


V 


m 


(3.1) 
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Control is accomplished by means of varying the magnitude and direction 
of spacecraft thrust, 



Thrust magnitude T is constrained consistent with rated maximum thrust 
and the number of working thrusters. Control of thrust direction [5 is constrained 
in such a way as to reflect the full 360 degree range of motion afforded by 
slewing the entire spacecraft. 

0<r<2x(92xl0^^A^) 

0< P<2n 

3. Dynamic Equations 

The dynamic equations used to model the reference mission trajectory are 
based on a heliocentric two body model [Refs. 8, 9] and can be derived using 
LaGrangian mechanics as follows. 

d (dL _ j 
dt ^dxj dx 

in which is the generalized non-conservative force vector representing 

spacecraft thrust. System kinetic and potential energy terms define the 
LaGrangean, 

T = -V-V = -(r^ + r^0^) 

2 2 

y_-GM_ // 
r r 

L = T-V = -ir^ + r^e^) + ^ 

2 r 

Spacecraft acceleration (thrust divided by mass) is described in polar 
coordinates, 
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T 

a = — 
m 


a = acos /36 + asm ( 5 r 


ttg = -COS P 


T 

=—sin (3 
m 


Generalized specific non-conservative force and torque vectors, in the 
radial and angular direction respectively, can then be defined as. 


= asin P 


Qg = ra cos P 


and applied to derive the dynamic equations such that. 


in which. 


such that. 


dtydr ] dr 


d ( dL 


dt V dO 


-- = 2 
d 0 


d f dL 


dt V dr 


—\^] = 2rW + r^0 
dtydO J 




v = rd 


v = rd 
v = rd + rO 


r - rO^ = asm P 
r 


IrrO + r ^6 = ra cos P 


■ V ju . „ 

u -+ -^ = asinp 

r r 


■ o 

V -I-= a cos p 

r 


■ ^ ■ o 

u = --^ + asinp 

r r 


n 

v = a cos p - 

r 
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which results in dynamic equations that fully describe the two body heliocentric 
model, 



v" // r . _ 

-—Sin P 

r r m 

-uv T 

-^—cos p 

r m 


( 0 . 1 ) 


where is thruster exhaust velocity and p is the solar gravitational constant. 

The gravitational effects of planets and asteroids can be neglected for the 
purposes of this simplified model without substantially affecting the outcome of 
the results. [Ref. 9] 

B. MODELING THE EVENTS 

To generate an optimal trajectory, it is first necessary to specify initial and 
final events, denoted by e, in terms of the five variables comprising the 
spacecraft state vector at any given point. Wherever this can not be done 
explicitly, an orbital element or other quantity must be calculated from the state 
variables. This element or quantity can then be constrained to correspond to the 
requirements of a particular event. Physically, these events represent launch 
from Earth, the landing on 1989ML departure from 1989ML and 

subsequent return to Earth ,e /. 
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It should be noted that the bold letter C denotes an event condition 
including all of its corresponding elements as shown above. Whereas the 
ordinary e simply denotes orbital eccentricity at a particular event. r,0,u,v and 
m are state variables and are subscripted above to denote their value at a 
particular event. 

C3 is traditionally taken to be the major performance parameter agreed on 
between the launch vehicle system and the planetary spacecraft [Ref. 10]. The 
injected mass or payload mass for a given launch vehicle and a given trajectory 
is always a function of C3. In other words, the spacecraft wet mass is dictated by 
C3. In this case, C3 equals the square of the hyperbolic excess velocity of the 
spacecraft with respect to Earth at the initial event, and constitutes the launch 
energy required to get to 1989ML using the given trajectory. It is therefore 
necessary to constrain C3 within the performance parameters of the booster as 
part of the initial event constraints. 

signifies final hyperbolic excess velocity with respect to Earth at the 

final event, a represents the semi-major axis of the spacecraft orbit at a 
particular event, at is the argument of periapsis (measured from the Ascending 
Node), d) is the longitude of periapsis (measured from the Vernal Equinox), t is 
the time, and M is the mean anomaly corresponding to a particular event. 

The time between the two interior events is spent on the asteroid surface 
and the spacecraft trajectory is assumed to be identical to that of the asteroid. 
Thus, it is not subject to optimization. A MATLAB application package called 
DIDO is used throughout this work to generate optimal trajectories based on the 
specified dynamics, events, and constraints. 

1. Initial Events 

For the purposes of this simplified analysis, the Earth’s orbit around the 
Sun is assumed to be perfectly circular with a radius of 1 Astronomical Unit or 1 
AU (approximately 1.495xl0®fcm). The initial spacecraft radial position is thus 
specified at 1 AU, but initial angular displacement is completely unconstrained 
leaving DIDO to choose the optimum angular position, and thus the 
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corresponding calendar date, from which to launch in order to minimize fuel 
consumption consistent with mission constraints. Spacecraft wet mass can be 
specified at the beginning of any scenario. This mass is then used as the initial 
value for the corresponding state variable. 

Maximum C3 can also be specified at the beginning of each scenario. 
DIDO will choose an initial spacecraft velocity vector so as to minimize ion 
thruster fuel consumption within the constraints of the maximum C3 specified. 
Initial event constraints are summarized below. 


0 < < 2;r 


- ^0 - 

“V^^max - ^0 - 
0 ^ Caused ^ 


In order to constrain C3 within the maximum C3 specified at the beginning 
of the simulation, it is first necessary to define C3 in terms of the spacecraft state 
variables. The magnitude of the initial spacecraft velocity (in scaled units or 
km/sec) is given by equation 3.5. 

K=V“o+Vo (3.5) 


Initial flight path angle is given by equation 3.6. 

fpa^^ = arctan(Mo / v„) (3.6) 


These terms are employed in equation 3.7 [Ref. 10] through which the 
spacecraft state vector can be constrained consistent with a specified maximum 
C3. 


C3 = KL + cos(^a„) (3.7) 

2. Interior Events 

Between the time that the spacecraft touches down on 1989ML and the 
time it lifts off, the spacecraft state variables representing position and velocity 
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will change consistent with the orbit of 1989ML Since the asteroid’s semi-major 
axis, eccentricity, and longitude of periapsis remain constant throughout the orbit, 
these quantities are calculated in terms of the spacecraft state variables and then 
constrained to equal the corresponding values for 1989ML This effectively 
collocates the spacecraft with 1989ML for the surface time specified at the 
beginning of a simulation. Since a primary motivation for this thesis is to 
determine whether or not changes in asteroid surface time can be used to reap 
significant savings in fuel mass, it is desirable to run different mission scenarios 
in which various asteroid surface times can be imposed. 

For a circular orbit, angular displacement is directly proportional to 
elapsed time. For elliptical orbits, mean anomaly is directly proportional to 
elapsed time. Thus, changes in mean anomaly are used here to solve for and 
constrain asteroid surface time within a specified window [Ref. 18]. 

To find mean anomaly in terms of the state variables using Kepler’s 
Equation, it is first necessary to calculate the semi-major axis, and eccentricity 
vector for the event in question. After converting spacecraft position and velocity 
into cartesion coordinates, 

rcosO MCOS^-vsin6* 

R= rsinO V= Msin^-vcos6* 

0 J L 0 

Specific mechanical energy is found in terms of spacecraft orbital radius 
[Ref. 7]. 

I -|2 

^ u 

f = LL_-fi (3,8) 

2 r 

Semi-major axis is then calculated [Ref. 7]. 



(3.9) 


Along with the eccentricity vector [Ref. 7]. 
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(3.10) 



Eccentric anomaly is then found using the arctangent of the following two 
expressions in which u is the true anomaly [Ref. 7]. 


sinE = 


siny 



1 + cos u 


(3.11) 


^ ae + rcoso 

cos E =- 

a 


(3.12) 


which is all that is needed to solve for mean anomaly at asteroid arrival and 
asteroid departure using Kepler’s Equation (Equation 3.14). Mean motion while 
on the asteroid is found using the semi-major axis and the solar gravitation 
constant in Equation 3.13. 


«= f4 (3.13) 

V a 

M=E-esmE (3.14) 

A corresponding time can thus be calculated for arrival or departure. 

t- — 
n 

Longitude of periapsis can be found using Equation 3.15 [Ref. 7], 

cos(®) = 4—^ (3.15) 

l\e\ 

taking advantage of the small angle approximation, Equation 3.15 can be 
constrained such that, 

a) = Q. + a) 

where the argument of periapsis of 1989ML is defined as ® = 183.28° and the 
longitude of the ascending node for 1989ML is defined as Q = 104.42°. 
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Finally, spacecraft mass at rendezvous is set equal to the spacecraft mass 
at departure since no thrusting should occur while the spacecraft is on the 
surface. 


3. Final Events 

Final spacecraft radius at Earth return must equal 1 AU since Earth orbit 
around the Sun is assumed to be perfectly circular. Once again, DIDO is free to 
pick final angular displacement, but subject to the constraint imposed by Earth 
asteroid phasing requirements. Thus, the angular position of Earth is moved 
forward from the initial launch position of the spacecraft in direct proportion to 
Earth mean motion and total flight time. The final spacecraft angular position is 
then constrained to this value such that. 

Unlike the rendezvous with 1989ML, in which spacecraft position and 
velocity are matched to the asteroid. Earth return can be characterized as an 
intercept in which some hyperbolic excess velocity (typically no more than 5 
kilometers per second) is permitted. This value can be specified at the beginning 
of each simulation. DIDO allocates the specified magnitude among the 
transverse and radial components of the final state in an optimal fashion; shaping 
the trajectory to minimize ion thruster fuel consumption. 

C. SIMPLIFIED MODEL OF SPACECRAFT AND THRUSTER 

The spacecraft model used in the simulations featured here assumes a 
fully fueled or “wet” mass of 1222 kg and mounts three NSTAR thrusters; two of 
which are operated at any given time. 
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NSTAR 

Throttle 

Level 

Mission 

Throttle 

Level 

PPU Input 
Power 
(kW) 

Engine 
Input Power 
(kW) 

Measured 

Thrust 

(mN) 

Main Flow 
Rate 
(seem) 

Cathode 
Flow Rate 
(seem) 

Neutralizer 
Flow Rate 
(seem) 

Speelfle 

Impulse 

(s) 

Total 

Effieieney 

12 

85 

1.99 

1.86 

75.34 

19.99 

") c ; 1 

2.82 

3035 

0.602 

1 1 

83 

1.94 

1.82 

72.55 

18.63 

ri 

2.67 

.3125 

0.610 


83 

1.96 

1.83 

72.63 


2.75 

2.67 


0.609 

10 

77 

1.84 

1.72 

69.54 


2.75 

2.67 


0.594 

10 

76 


1.70 


17 . 31 

2.58 

2.51 



10 

75 


1.68 


17.33 

2.58 

2.51 


0.601 1 

10 

74 

1.77 

1.66 

66.11 


2.59 

2.51 

MtTifB 


10 

73 

1.75 

1.65 

65.64 


2.59 

2.51 



10 

72 

1.73 

1.63 

65.15 

17 . 31 

2.59 

2.51 



9 

69 

1.67 

1.57 

62.27 

16.08 

2.50 

2 . 4.3 

3070 

0.597 

6 

48 

1.29 

1.22 

47.43 

1 1.42 

2.50 

2.42 

3006 

0.573 

6 

48 

1.29 

1 O') 

47.39 

1 1.44 

2.49 

2.42 

3004 

0.571 

3 

27 

0.89 

0.84 

31.70 

6.93 

2.50 

2.43 

2770 

0 . 511 

0 

6 

0.50 

0.48 

20.77 

6.05 

2.50 

2.43 

1961 

0.418 


Table 1. Flight Engine Performance Measured in Space [From Ref. 12] 


To model the performance of the NSTAR ion thruster, values for 
Measured Thrust and corresponding Specific Impulse are selected or 
interpolated from Table 1. While the actual NSTAR thruster is capable of 
automatically switching between dozens of Mission Throttle Levels to 
accommodate the fluctuating availability of electrical power, the simplified model 
used here neglects to account for this feature; nor does it account for variations 
in solar flux levels as a function of distance from the Sun. Therefore, Specific 
Impulse and Maximum Thrust must be specified at the beginning of each 
simulation. During the simulation, spacecraft thrust ranges all the way from zero 
to the specified maximum thrust constraint while specific impulse remains 
constant at all times. 

D. VALIDATING THE MODEL 

1. Optimality 

After a mission trajectory is generated, optimality can be verified using the 
basic principles of optimal control theory. The theory states that in order for a 
trajectory to be optimal, the Flamiltonian must equal zero at all points along the 
trajectory. Furthermore, the derivative of the LaGrangean of the Flamiltonian with 
respect to the individual controls should also be zero while fulfilling the 
requirements of the Complimentarity Condition [Ref. 8, 13]. 
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The Hamiltonian of the dynamics described by equation 3.3 is defined 
such that: 


H (^A,x,u,t) = F (^x,u,t) +f (3.16) 


in which A,x,u,t represent the costate, state, control, and time respectively. 
Here, optimality is defined only in terms of an endpoint cost - final spacecraft 
mass. Therefore, the first term in the right hand side of the equation equals zero. 
The second term is expanded to show all state and costate variables below. 


,v v JU T . uv T n n T 

H =A^u+\- + \ -+ T„—sin^-T, —+ T,—cos^-4 — 

r r r m r m 


(3.17) 


Applying the Karush-Kuhn-Tucker (KKT) condition in which. 


dH 

du 


_ dH 
du 


dH 

du 


ju = 0 


where the LaGrangian of the Hamiltonian can be defined as; 

H (^A,x,u,t,ju) = H (^A,x,u,t) + f/h{u,t) (3.18) 

in which ju is the KKT multiplier and h is the constraint function [Ref. 13].and 
taking the derivative with respect to each control variable yields, 

r)H 11 1 

^ = 4—+ cos^-4 —+ (3.19) 

Cl m m V 

e 


r)M T T 

= + (3.20) 

dp m m 

If in fact the generated trajectory is optimal, the Minimum Principle 
requires that the results of equations 3.17, 3.19, and 3.20 be close to zero at any 
point in the trajectory. While this does not assure optimality, a result substantially 
greater than zero usually proves that the solution is not optimal. 

Further validation is provided by exploiting the Complementarity Condition 
[Ref. 13]. The KKT multiplier for each of the two control variables on the left in 
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Equations 3.21 and 3.22 must behave according to the control inequalities on the 
right for a potential optimal trajectory to exist. 

^o^r = r_ 

=o^7;„<r<r_ (3.21) 

>0^T = T 

V./ ' 

<0^ P = 271 

ju^<=0^0</3<27r (3.22) 

[ >0^/3 = 0 

2. Feasibility 

Generated trajectories are checked for feasibility by taking the control 
profile created by DIDO and propagating it through the equations of motion 
shown in equation 3.3 using the ODE 45 propagator in MATLAB and linear 
interpolation of the DIDO controls. The propagator trajectory is then 
superimposed on the optimizer trajectory. A perfect match indicates feasibility. 

3. A Word about Scaling 

Ouantities calculated in the simulations featured here range from the very 
large (interplanetary distances) to the very small (spacecraft thrust). While 
employing standard units of measure may be convenient for the reader, such 
units can be computationally inefficient or even inaccurate. Thus, it becomes 
necessary to define an alternative method of scaling for fundamental units such 
as time, distance, and mass. Done properly, this can balance out the calculations 
in such a way as to ensure well behaved solutions and computational accuracy. 
These fundamental units are then used to generate derived units such as 
velocity, acceleration, and thrust. For example, it is convenient to define the unit 
of mass in terms of the initial spacecraft mass, distance in terms the distance 
between the Earth and the Sun, and time in terms of a constant force unit. 

Mass Unit = 1222 kg = 1 MU 
Distance Unit = 1.49598x 10^A:m = 1 DU 
Force Unit = 0.0145^^-fcm/sec^ = 1 FU 
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Time Unit = 

V FU 

Velocity Unit = DU/TU 

It is occasionally desirable to readjust or rebalance the relationship 
between all of these units. This is particularly true when the predefined values 
like rated maximum thrust or spacecraft mass change significantly. Inadequate 
scaling manifests itself in a number of ways such as a large Hamiltonian, 
inconsistent propagation, or even through an altogether infeasible solution. The 
author found it convenient to address this issue by adjusting the value of FU 
between 0.0145 (for the Hohmann Transfer example) and 0.0148 (for all other 
examples) until the propagation inconsistencies subsided. 

In many of the examples that follow, state variables are plotted in both 
scaled and unsealed (metric) units for the reader’s convenience. If no units are 
specified, scaled units can be assumed. It should be noted that all the scaled 
plots that follow display horizontal axes which are numbered in the Time Unit 
described above. Each Time Unit is approximately equal to 41 days. 
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IV. SIMULATING THE ASTEROID SAMPLE RETURN PROBLEM 


A. ONE-WAY ORBIT TRANSFERS 

Simulation of the sample return mission is an incremental process in 
which increasingly complex trajectories are modeled, generated, and validated 
for optimality and feasibility at each step. Once simple trajectories are found to 
be valid, additional complexities can be encoded into the model that allow the 
simulation to account for such factors as C3, hyperbolic excess velocity, or 
eccentric orbits. Once complex one-way trajectories are validated, the code used 
to generate them is employed to create increasingly complex two way models 
that ultimately incorporate phasing, as well as rendezvous and capture. Where 
possible, additional validation is sought by way of comparison with trajectories 
generated by other researchers. 

1. Simple Circular Transfer 

A simple Hohmann Transfer was chosen to validate the dynamics and 
ensure the soundness of the code and scaling. Using a known spacecraft model 
and a known optimal trajectory with known constraints, an attempt was made to 
duplicate the results obtained by Lieutenant Scott B. Josslelyn [Ref. 9] and 
thereby confirm the code’s accuracy. The model assumes a 659.3 kilogram 
spacecraft powered by six NSTAR thrusters with a total maximum rated thrust of 
0.55 N and a specific impulse of 3280 seconds. The trajectory is constrained to 
take the spacecraft from an assumed Earth circular orbit at 1 AU to an assumed 
Mars circular orbit at 1.52 AU while traveling through no more than 3.14 radians 
of angular displacement. 

Where Jossleyn obtained a final spacecraft mass of 540.3 kilograms, the 
sample return mission simulator code produced a final mass of 540.7 kilograms. 
Similarly, where Jossleyn obtained an elapsed time of 253.4 days, the asteroid 
code generated an elapsed time of 252.8 days. The trajectory shown in Figure 7 
is virtually an exact replica of the one obtained by Jossleyn and is similar in 
appearance to a Hohmann Transfer. 
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In addition, the control and state variables plotted in Figure 8 and Figure 9 
exhibit smooth, economical, non-erratic behavior with respect to time. Thrust 
magnitude exhibits a profile characteristic of “bang-bang control” as one would 
expect in a Flohmann Transfer. These are all traits highly suggestive of a solution 
that is feasible if not optimal. 
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Figure 8. Hohmann Transfer - State and Thrust Angle (Scaled Units) 
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Figure 9. Flohmann Transfer - State and Thrust Magnitude (Unsealed Units) 
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Optimality is apparent however, in the plots of the Hamiltonian, KKT 
multipliers, and switching functions with respect to time shown in Figure 10. 



Figure 10. Hohmann Transfer - Hamiltonian and KKT Multipliers 

(Scaled Units) 


That the KKT multiplier for thrust varies somewhat from zero may be due 
to a problem within DIDO requiring further investigation [Ref. 15]. Non-optimal 
solutions frequently have KKT thrust values many orders of magnitude higher 
than the one shown above. This plot suggests that the solution is at least very 
nearly optimal. 

But taken in conjunction with Figure 11, which shows the Complimentarity 
Condition being met throughout the trajectory. Figure 10 is highly indicative of an 
optimal solution. 
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Figure 11. Hohmann Transfer - Switching Functions (Scaled Units) 


Finally, feasibility can be inferred based on the appearance of Figure 12; 
where state trajectory, shown in red, is superimposed on propagated trajectory, 
shown in blue. The two are virtually indistinguishable. 

Thus, the trajectory is in all likelihood optimal and feasible; exhibiting 
behavior perfectly consistent with the necessary conditions for optimality. 
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Figure 12. Hohmann Transfer Feasibility Comparison (Scaled Units) 

2. Dealing with the Limitations of the Simplified Model 

The previous trajectory was generated for a 660 kilogram spacecraft with 
six thrusters operating at full power. For the purposes of this work, the sample 
return spacecraft is assumed to be a 1222 kilogram spacecraft with two thrusters 
operating at less than full power much of the time. 

As previously discussed, a real NSTAR ion thruster features a range of 
mission throttle levels; some of which are shown in Table 1. Each level has a 
particular thrust setting and a corresponding specific impulse. The choice of what 
throttle level to employ is dictated not only by trajectory requirements, but also by 
available power levels. These are in turn a function of available solar flux, battery 
and solar array condition, electrical loads, and many other factors. 

None of this is adequately accounted for in the simulations that follow. The 
thruster model used here is simplistic and artificial. Specific impulse is fixed 
throughout the duration of the mission simulation and thrust is allowed to sweep 
between zero and a designated maximum value irrespective of real-world power 
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constraints. Moreover, the sweep takes place in a continuous fashion, not in the 
discrete fashion exhibited by NSTAR. All of this gives rise to the question of how 
best to account for these artificialities within the limitations of the model. 

It is unrealistic to suppose that a real spacecraft will employ maximum 
thrust or be able to operate at maximum specific impulse at all times due to 
constraints on available power. It is equally unrealistic to suppose that a real 
spacecraft will always employ zero or negligible thrust over the course of a 
voyage since some amount of thrusting is necessary to arrive at the destination. 
Nevertheless, the inherent limitations of the model force one to choose from 
between the two extremes. 

Thus, it was initially decided to split the difference and assume a 
maximum thrust about halfway between zero and the maximum thrust value 
shown at the top of Table 1. The associated specific impulse was obtained by 
interpolating between the specific impulse values corresponding to throttle levels 
three and six. While the numbers may not correspond precisely to an actual 
Mission Throttle Level, they are close enough for the purposes of this analysis. 

Subsequent experimentation revealed that increasing the thrust by about 
one mN to 38.21 mN and changing the specific impulse proportionally to 2867.2 
seconds for a given set of mission constraints (including 90 days on the asteroid) 
resulted in an overall mission fuel expenditure similar to that obtained by the JPL 
software package known as SEPTOP. This is not to suggest that the results of 
the simulations discussed in this thesis exhibit anything approaching the 
accuracy or fidelity of SEPTOP. The highly simplified two body model used here 
doesn’t even begin to compare. 

Picking a maximum thrust / specific impulse combination that falls about 
half way between the two demonstrated performance extremes of the real 
thruster and that also happens to generate a result roughly consistent with a 
much higher fidelity simulation might seem crude. But it is also very conservative 
and should go some way toward compensating for the artificialities mentioned 
above. 
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3. Earth to Asteroid Transfer 

In addition to using thrust and specific impulse values that fall within 
NSTAR performance parameters, the Earth to Asteroid transfer improves upon 
the Hohmann Transfer by modeling the orbit of 1989ML as an ellipse and by 
providing an allowance for C3 to be imparted to the spacecraft by the booster. 
For this simulation, C3 was set not to exceed 2.67 W/sec^ so as to emulate the 
performance of the Delta II launch vehicle [Ref. 16] which can actually launch 
1222 kg to a C3 of 2.67 W/sec^ The transfer takes 1.56 years and consumes 
67.7 kilograms of Xenon. The resultant optimal trajectory and thrust vector profile 
is shown below. 



It is worth noting that the thrust vectors have a comparatively small radial 
component. C3 has provided much of the energy needed to accomplish the 
transfer so the thrusters do not need to contribute that much more. In later 
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simulations it will be shown that all other things being equal, thrust direction 
becomes more radial as thrust magnitude is constrained. 



Figure 14. Earth to Asteroid Transfer - State and Thrust Angle (Scaled Units) 
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Figure 15. Earth to Asteroid Transfer - State and Thrust Magnitude 

(Unsealed Units) 


Once again, state variables plotted in Figure 14 and Figure 15 appear to 
exhibit smooth behavior and “bang-bang” control. The Flamiltonian and KKT 
multipliers are very close to zero. Minor deviations from zero on the part KKT 
Thrust are due once again to a preexisting problem within DIDO requiring further 
investigation. A quick glance at Figure 17 verifies that complimentarity has been 
satisfied and a comparison of the propagated trajectory with the proposed 
optimal trajectory in Figure 18 shows that the trajectory is feasible. 
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Figure 16. Earth to Asteroid Transfer - Hamiltonian and KKT Multipliers 

(Scaled Units) 
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Figure 17. Earth to Asteroid Transfer - Switching Functions (Scaled Units) 
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Figure 18. Earth to Asteroid Transfer - Feasibility Comparison (Scaled Units) 

4. Asteroid to Earth Transfer 

The second half of the asteroid sample return mission shown below 
entails Earth intercept from elliptical orbit. Flyperbolic excess velocity at Earth is 
now incorporated into the final event and set not to exceed 3.94 km / sec. The 
transfer takes 0.84 years to accomplish and consumes 27.35 kilograms of 
Xenon. Modeling Earth arrival as an intercept as opposed to a rendezvous 
equates to a considerable savings in fuel since less thrusting is necessary. 

As in the previous two trajectories, evidence of optimality is found by virtue 
of the Minimum Principle and Complimentarity. Results are plotted in Figure 22 
and Figure 23. Feasibility is clearly in evidence in Figure 24. 
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Figure 19. Asteroid to Earth Transfer - Trajectory (Unsealed Units) 



Figure 20. Asteroid to Earth Transfer - State and Thrust Angle (Scaled Units) 
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Figure 21. Asteroid to Earth Transfer - State and Thrust Magnitude 

(Unsealed Units) 
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Figure 22. Asteroid to Earth Transfer - Hamiltonian and KKT Multipliers 

(Scaled Units) 
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Figure 23. Asteroid to Earth - Transfer Switching Functions (Scaled Units) 
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Figure 24. Asteroid to Earth Transfer - Feasibility Comparison (Scaled Units) 
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B. TWO WAY ORBIT TRANSFERS 

1. Unconstrained Surface Time with No Phasing 

Having successfully modeled and validated trajectories from Earth to the 
asteroid and from the asteroid to the Earth, the next logical step is to incorporate 
everything used thus far into a single round trip simulation. The resulting optimal 
trajectory allows the spacecraft to stay in the orbit of 1989ML for 355 days. It 
cannot actually be claimed that the spacecraft has rendezvoused with the 
asteroid and returned to Earth since phase differences between Earth and 
1989ML have yet to be taken into account. That step comes later. For this 
trajectory, shown in Figure 25, DIDO has elected to use only 2.6 of the allowed 
3.94 km / sec in hyperbolic excess velocity. Total fuel consumption is 93 
kilograms and total mission time is 3.98 years. It is also worth noting that DIDO 
has elected to place the points of launch and return to Earth near asteroid 
perihelion while placing the asteroid arrival and departure points in the vicinity of 
asteroid aphelion. There is almost a symmetry in the trajectory and thrust profile 
that mirrors the symmetry of the 1989ML orbit itself. The outbound portion of the 
trajectory is denoted in blue. The inbound or return trajectory is denoted in red. 



Figure 25. Unconstrained Surface Time with no Phasing - Trajectory 

(Unsealed Units) 
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The discontinuities appearing in the state variable profiles of Figure 26 
and Figure 27 mark the point at which trajectory optimization stopped and then 
started up again. These were described earlier as the left and right interior events 
between which the spacecraft is said to be co-orbiting with 1989ML such that no 
thrusting occurs and no trajectory optimization is necessary. 

Figures 28 and Figure 29 show the propagated trajectory superimposed 
on the optimal trajectory for the outbound and inbound transfers respectively. 
Separating the two halves of the trajectory makes it easier to see that the 
trajectories appear to be virtually identical. This once again serves to verify 
feasibility. 



Figure 26. Unconstrained Surface Time with no Phasing - State and Thrust 

Direction (Scaled Units) 
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Figure 27. Unconstrained Surface Time with no Phasing - State and Thrust 

Magnitude (Unsealed Units) 



Figure 28. Unconstrained Surface Time Outbound - Feasibility Comparison 

(Scaled Units) 
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Figure 29. Unconstrained Surface Time Inbound - Feasibility Comparison 

(Scaled Units) 


2. Unconstrained Surface Time with Phasing 

The next step involves using mean anomaly as discussed previously to 
impose the condition that the spacecraft return to Earth orbit at the point where 
Earth is actually located at the end of the mission. The first thing to notice is that 
fuel consumption increases dramatically to 156 kilograms of Xenon. Unlike the 
previous case in which phasing was not accounted for; this time the spacecraft 
can actually be said to have landed on 1989ML and returned to Earth. Flowever 
the phasing constraint permits a stay of only 28 days. A longer stay would require 
more fuel and be less efficient. As in the previous example, discontinuities in the 
state variable plots show that no optimization is taking place while the spacecraft 
is on the asteroid surface. Once again, the outbound and inbound trajectory 
propagations confirm the feasibility of the solution as shown in Figure 33 and 
Figure 34. 


43 


























Figure 30. Phasing - Trajectory (Unsealed Units) 



Figure 31. Phasing - State and Thrust Direction (Scaled Units) 
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Figure 32. Phasing - State and Thrust Magnitude (Unsealed Units) 
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Figure 33. Phasing - Outbound Feasibility Comparison (Scaled Units) 
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Figure 34. Phasing - Inbound Feasibility Comparison (Scaled Units) 

3. Constrained Surface Time with Phasing 

The final step in modeling the full sample return is to enable the mission to 
be simulated with asteroid surface time constrained to a specific value or range 
of values. As discussed earlier, this is done by making use of the left and right 
mean anomalies calculated for each of the two interior events. 

Figure 35 shows this comprehensive model trajectory for a 90 day 
asteroid visit. Total mission time is 3.17 years. C3 is 2.67 W/sec^ and 
hyperbolic excess upon Earth return is 3.2 km / sec. Total mission fuel 
consumption is 180.13 kilograms of Xenon. 

Figure 36 and Figure 37 show the associated scaled and unsealed state 
and control variables while Figure 38 and Figure 39 illustrate that the trajectory is 
in fact feasible. 
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Figure 36. Comprehensive Model - State and Thrust Magnitude 

(Unsealed Units) 
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Figure 37. Comprehensive Model - State and Thrust Direction (Scaled Units) 
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Figure 38. Comprehensive Model Outbound - Feasibility Comparison 

(Scaled Units) 
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Figure 39. Comprehensive Model Inbound - Feasibility Comparison 

(Scaled Units) 


4. Comparison with SEPTOP Resuit 

Comparing the results of a high fidelity, n-body, three dimensional model 
like SEPTOP with the results of the highly simplified 2-body, two dimensional 
model can be misleading. Nevertheless, there are interesting similarities between 
the SEPTOP result shown in Figure 40 and the Comprehensive Model result 
shown in Figure 35. Asteroid surface time is the same. The arrival and departure 
points between orbits occur at roughly the same positions. Total flight time is 
roughly the same. C3 is the same as are wet mass and propellant mass. There 
are lesser similarities in hyperbolic excess velocity and outbound thrust profile. 

Fixed values for maximum thrust and associated specific impulse were 
extrapolated from Table 1 in such a way as to result in a fuel expenditure 
consistent with Figure 40 while keeping near the middle of the thrust range. 
Nevertheless, the fact that all the other similarities mentioned above fell in to 
place at the same time suggests that even the highly simplified model might be 
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sophisticated enough to prove useful in formulating fast rough estimates of how 
changes in asteroid surface time can affect required fuel load. 

1989 ML 

Delta 2925-9.5 / SEP (6.3 kW) 



1989 ML rendezvous 
12-2-10 Tw=90d 


30 day tics on spacecraft path 


Mars 


Earth Return 7-7-12 
Vhp= 3.94 km/s 


3.15 Year Flight Time 
0^= 2.67 km^/s^ 
DL/^ -27.0 deg 
Mwet= 1 222 kg 
MP=180 kg 
MDRY= 1042 kg 
AV= 4.829 km/s 
2 operating NSTAR 
Thrusters 

6.3 kW EOL Solar Array 
494 W Spacecraft Power 
0% LV Contingency 
95% SEP Duty Cycle 


Sauer 2-1804 
1989n1c 


Figure 40. SEPTOP Trajectory Diagram (Courtesy of JPL) 


C. VARYING THE SURFACE TIME AND EFFECT ON FUEL 

A maximum thrust of 35.14 mN and a corresponding interpolated specific 
impulse of 2821.6 seconds was used to generate the top line in Figure 41. The 
combination constitutes the lowest throttle level that can be employed with this 
code to yield feasible results over a significant range of asteroid surface times. 
Similarly, a maximum thrust of 75.34 mN and a corresponding specific impulse of 
3035 seconds (Throttle Level 12 from Table 1) was used to generate the bottom 
line. The “R squared” value denotes the square of the correlation coefficient 
between the points plotted and the equation of the straight line drawn through 
them. Perfect correlation in the linear regression of the top two lines would result 
in a value of one. This won’t be true of the bottom line since the slope is very 
close to zero. But the correlation is evident at a glance. 
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Thus, the top and bottom lines, derived by testing the demonstrated 
operating extrema of the DS 1 NSTAR thruster as described in Table 1, 
constitute best and worst case relationships. At worst, staying an extra day on 
1989ML will cost no more than 0.64 kilograms of fuel. At best, it will cost next to 
nothing. Most likely, it will cost something in between. The middle line is derived 
from the maximum thrust / specific impulse combination that matched the 
SEPTOP results. It suggests that additional time spent on the asteroid can be 
purchased at a cost of 0.5 kilograms of fuel per day. 

Attempts to extend these lines further to the right resulted in infeasibilities. 
This suggests that at some point, no amount of fuel will get the spacecraft back 
to Earth on time. However, High specific impulse and high thrust certainly 
increases one’s options and flexibility. 



D. A TENDENCY TOWARD RADIAL THRUSTING 

In the course of research, it was frequently observed that optimal thrust 
angle trajectories for spacecraft that needed to thrust a great deal in order to 

complete a transfer tended to become more radial as the spacecraft approached 
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aphelion. At first, this was attributed to a plotting error since radial thrusting does 
not strike this author as being particularly efficient much less optimal. 

The phenomena can be observed in the following trajectory plots. The 
same code was used to generate Figures 42, 43, and 44. They show 
comprehensive model trajectories for asteroid sample return missions of 106, 
110, and 100 day surface times respectively. The only difference between them 
is that maximum available thrust and specific impulse is progressively 
constrained from Figure 42 onward. Each trajectory was found to be both optimal 
and feasible. 
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Figure 43. 47.43 N at 3006 sec (110 day surface time) 



Figure 44. 35.14 N at 2821.6 sec (100 day surface time) 
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V. CONCLUSIONS AND FUTURE WORK 


It would appear, based on Figure 41, that fuel consumption over the 
course of a sample return mission optimal trajectory is a linear function of 
asteroid surface time; at least within the limitations of the thruster model 
employed here. For the assumed case, time spent on the asteroid is purchased 
at a cost of about 0.5 kilograms of Xenon per day. While the assumed case 
maybe overly conservative, similarities between it and the results shown in 
Figure 40 suggest that it might be a good starting point from which to extract 
quick rough estimates of how changes in asteroid surface time might affect fuel 
load. Of course, the higher the available thrust and specific impulse, the more 
scheduling flexibility one has in exploring the surface of 1989ML. If the spacecraft 
can be designed to sustain thruster performance consistent with NSTAR Throttle 
Level 12 (The lower line on Figure 41) over the course of the entire mission, then 
asteroid surface time can range from forty to one hundred and sixty days without 
any appreciable cost in additional fuel mass. 

Furthermore, as shown in Figures 42 through 44, the optimal thrust angle 
trajectory for the sample return spacecraft becomes increasingly radial as 
thruster performance becomes increasingly constrained. 

Flopefully, these insights will be of use to mission planners at JPL and 
elsewhere. As ion thrusters continue to proliferate, questions of optimal control 
will become more important and there are many aspects of the topic which need 
to be explored. Among these are; 

Three Dimensional Orbit Transfer Models. The author has put 
considerable effort into formulating a MATLAB code to generate optimal low 
thrust trajectories in three dimensions using spherical coordinates; with little 
success thus far. One of the first steps in the validation process involved 
generating coplanar Flohmann transfers at various angles of inclination to ensure 
that all the results were identical. Lovely, well behaved, three dimensional 
coplanar transfers were generated only when available fuel mass was artificially 
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restricted to that quantity known to be needed for an optimal (Hohmann) transfer. 
Left to it’s own devices by removal of the fuel mass constraint, the code 
generated a hideously incontinent thrust profile in which the spacecraft weaved in 
and out of plane in a manner than was clearly less than optimal. While the 
spherical dynamic equations are subject to singularities under certain 
circumstances, no definitive cause for the difficulties could be established at the 
time of this writing. 

A Realistic NSTAR Thruster Model. The NSTAR thruster is likely fly on a 
number of different spacecraft for a long time to come. Some way of realistically 
coding and modeling the NSTAR throttle settings and accounting for variations in 
solar flux as a function of solar distance would go a long way toward increasing 
the fidelity of this and other models. 
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